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Abstract 

In the study of ocean wave impact on structures, one often uses Froude scaling since 
the dominant force is gravity. However the presence of trapped or entrained air in 
the water can significantly modify wave impacts. When air is entrained in water in 
the form of small bubbles, the acoustic properties in the water change dramatically. 
While some work has been done to study small-amplitude disturbances in such 
mixtures, little work has been done on large disturbances in air- water mixtures. We 
propose a basic two-fluid model in which both fluids share the same velocities and 
analyze some of its properties. It is shown that this model can successfully mimic 
water wave impacts on coastal structures. The governing equations are discretized 
by a second-order flnite volume method. Numerical results are presented for two 
examples: the dam break problem and the drop test problem. The results suggest 
that this basic model can be used to study violent aerated flows, especially by 
providing fast qualitative estimates. 

Key words: free-surface flow, wave impact, two-phase flow, compressible flow, 
finite volumes 



1 Introduction 



One of the challenges in Computational Fluid Dynamics (CFD) is to deter- 
mine efforts exerted by waves on structures, especially coastal structures. The 
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flows associated with wave impact can be quite complicated. In particular, 
wave breaking can lead to flows that cannot be described by models like e.g. 
the free-surface Euler or Navier-Stokes equations. In a free-surface model, the 
boundary between the gas (air) and the liquid (water) is a surface. The liq- 
uid flow is assumed to be incompressible, while the gas is represented by a 
medium, above the liquid, in which the pressure is constant (the atmospheric 
pressure in general). Such a description is known to be valid for calculating 
the propagation in the open sea of waves with moderate amplitude, which do 
not break. Clearly it is not satisfactory when waves either break or hit coastal 
structures like offshore platforms, jetties, piers, breakwaters, etc. 

Our goal here is to investigate a relatively simple two-fluid model that can 
handle breaking waves. It belongs to the family of averaged models, in the 
sense that even though the two fluids under consideration are not miscible, 
there exists a length scale e such that each averaging volume (of size e^) con- 
tains representative samples of each of the fluids. Once the averaging process is 
performed, it is assumed that the two fluids share, locally, the same pressure, 
temperature and velocity. Such models are called homogeneous models in the 
literature. They can be seen as limiting cases of more general two-fluid models 
where the fluids can have different temperatures and velocities [13]. Let us 
explain why it can be assumed here that both fluids share the same temper- 
atures and velocities. There are relaxation mechanisms that indeed tend to 
locally equalize these two quantities. Concerning temperatures, these are dif- 
fusion processes and provided no phenomenon is about to produce very strong 
gradients of temperature between the two fluids like e.g. a nuclear reaction in 
one of the two fluids, one can assume that the time scale on which diffusion 
acts is much smaller than the time scale on which the flow is averaged. Sim- 
ilarly, concerning the velocities, drag forces tend to locally equalize the two 
velocities. Deflne a time scale built on the mean convection velocity and a 
typical length scale. For flows in which the mean convection velocity is mod- 
erate, this time scale based on convection is much larger than the time scale 
on which velocities are equalized through turbulent drag forces. Hence, in the 
present model, the partial differential equations, which express conservation 
of mass (1 per fluid), balance of momentum and total energy, read as follows: 



where the superscripts ± are used to denote liquid and gas respectively. Hence 
a+ and a~ denote the volume fraction of liquid and gas, respectively, and 
satisfy the condition + a~ = 1. We denote by p^, tT, e respectively the 
density of each phase, the velocity, the pressure, the speciflc internal energy. 
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g is the acceleration due to gravity (in two space dimensions, g is equal to 
(0, —g))) P '•= c^^p^ + (y~p~ is the total density, E = e + is the specific 
total energy, H := E + p/p is the specific total enthalpy. In order to close 
the system, we assume that the pressure p is given as a function of three 
parameters, namely a = — a~ ^ p and e: 

p = V{a,p,e) . (5) 

We shall discuss in Section [2] how such a function V is determined once the two 
independent equations of state p = V^{p^^ e^) are known. Equations ([I])-([5j) 
form a closed system that we shall use to simulate aerated fiows. 

The main purpose of this paper is to promote a general point of view, which 
may be useful for various applications dealing with violent aerated fiows in 
ocean, offshore, coastal and arctic engineering. We do not consider here un- 
derwater explosions, where the word violent has a diflFerent meaning. The det- 
onation of an explosive charge underwater results in an initial high-velocity 
Shockwave through the water, in movement or displacement of the water itself 
and in the formation of a high-pressure bubble of high-temperature gas. This 
bubble expands rapidly until it either vents to the surface or until its internal 
pressure is exceeded by that of the water surrounding it [15]. What we do 
is to follow the approach first used, we believe, by the late Howell Peregrine 
and his collaborators [3lll7lll6j . The infiuence of the presence of air in wave 
impacts is a difiicult topic. While it is usually thought that the presence of 
air softens the impact pressures, recent results show that the cushioning eflFect 
due to aeration via the increased compressibility of the air-water mixture is 
not necessarily a dominant eflFect [1]. First of all, air may become trapped 
or entrained in the water in diflFerent ways, for example as a single bubble 
trapped against a wall, or as a column or cloud of small bubbles. In addition, 
it is not clear which quantity is the most appropriate to measure impacts. For 
example some researchers pay more attention to the pressure impulse than 
to pressure peaks. The pressure impulse is defined as the integral of pressure 
over the short duration of impact. A long time ago, Bagnold [Ij noticed that 
the maximum pressure and impact duration differed from one identical wave 
impact to the next, even in carefully controlled laboratory experiments, while 
the pressure impulse appears to be more repeatable. For sure, the simple one- 
fiuid models which are commonly used for examining the peak impacts are no 
longer appropriate in the presence of air. There are few studies dealing with 
two-fiuid models. An exception is the work by Peregrine and his collaborators. 
Wood et al. [21] used the pressure impulse approach to model a trapped air 
pocket. Peregrine & Thais [18] examined the effect of entrained air on a partic- 
ular kind of violent water wave impact by considering a filling fiow. Bullock et 
al. [5] found pressure reductions when comparing wave impact between fresh 
and salt water, due to the different properties of the bubbles in the two fiuids. 
Indeed the aeration levels are much higher in salt water than in fresh water. 
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Bredmose [2] recently performed numerical experiments on a two-fluid system 
which has similarities with the one we will use below. 

The novelty of the present paper is not the finite volume method used below 
but rather the modelling of two-fluid flows. Since the model described below 
does not involve the tracking nor the capture of a free surface, its integration is 
cheap from the computational point of view. We have chosen to report here on 
the case of inviscid flow. Should the viscosity effects become important, they 
can be taken into account via e.g. a fractional step method. In fact, when 
viscous effects are important, the flow is easier to capture from the numerical 
point of view. 

The paper is organized as follows. Section [2] provides an analytical study of 
the model. Section [3] deals with numerical simulations based on this model 
via a flnite volume method. Two examples are shown: the dam break problem 
and the drop test problem. Finally a conclusion ends the paper. 



2 Analytical study of the model 

2. 1 The extended equation of state 

It is shown in this section how to determine the function V{a^ p, e) in Eq. ([5j) 
once the two equations of state p = 7^^(p^, e^) are known. We call Eq. ([5]) an 
extended EOS, since 7^(— 1, p, e) = V~{p^ e) and p, e) = V^{p^ e), where 

p±=P±(p^e±), r± = T±(p±,e±), (6) 

are the EOS of each fluid, with the temperature of each phase. Although 
our approach is totally general, we will use the following prototypical example 
in this paper. Assume that the fluid denoted by the superscript — is an ideal 
gas: 

p- = (^-_l)p-e-, e- = CyT-, (7) 
while the fluid denoted by the superscript + obeys the stiffened gas law [GlfTT] : 

p+ + vr+ = (7+-l)p+e+, e+ = C+r+ + 4^, (8) 

where 7^, Cy, and tt"*" are constants. For example, pure water is well described 
in the vicinity of the normal conditions by taking 7"*" = 7 and tt"'" = 2.1 x 10^ 
Pa. 

Let us now return to the general case. In order to find the function there 
are three given quantities: « G [— 1, 1] , p > and e > . Then one solves for 
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the four unknowns p , e the foUowing system of four nonhnear equations: 



(l + tt)p+ + (l-tt)p-=2p, (9) 

(l + tt)p+e+ + (l-tt)p-e-=2pe, (10) 

P+(p+,e+)-7'-(p-,e-)=0, (11) 

T+(p+,e+)-T-(p-,e-) = 0. (12) 

For given values of the pressure p > and the temperature T > 0, we denote 
by 7^=^(p, T) and S^{p, T) the solutions (p=^, e^) to: 

7'±(p±,e±)=p, r±(p±,e±) = r, (13) 

and then: 

P = ^^^(^' + ^^"(P' , (14) 

pe = ^7e+(p, T) S^{p, T) + ^^7e-(p, T) £-{p, T) . (15) 

Finally the inversion of this system of equations leads to p = V{a^p^e) and 
T = T(a,p, e). 



Remark 1 T/ie system cind H^) is a differential and alge- 

braic equation, while the system ^\)-^ is a partial differential equation as it 
is the case for a system of single fluid equations. 

Concerning the prototypical case, the following generalization of ([71) is consid- 
ered: 

p- + 7r- = (7- - l)p-e-, e- = CyT' + . (16) 

7 P 

This generalization, which has the additional parameter tt", allows one to set 
the speed of sound to a certain value independently of 7", p~ and p~ . Using 
computer algebra to invert (fT^ and (fT5l) leads to the following expressions: 



7^(a,p,e) = (7(a) - l)pe-7r(a), (17) 

T(a.p.e)= ''^-'^"'°>:; 

pCy(a) 

where the five functions 7(a), 7r(Q;), Cv(a) and A^(q;) are defined by 
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' + ''+^, (19) 



7(tt) — 1 7+ — 1 7~ — 1 ' 
27r(«) _! + «_+, ^20) 

(21) 

7 ' — 1 7 — 1 

One can easily check that one recovers the equations of state for each fluid in 
the hmits a ^ ±1. Note that similar expressions can be found in Section 1.1 
of [12] where a two-dimensional, compressible, two-fluid mathematical model 
was used to compute numerically wave breaking. 




2.2 A hyperbolic system of conservation laws 



In this section, we assume that the system of equations is solved in M^, having 
in mind the numerical computations performed below. However the extension 
to 3D is straightforward. The system ([I])-(|4]) can be written as 

^ + V ■ ^(w) = 5(w) , (23) 

where 

w = := pui, pU2, pE) , (24) 

and, for every n G M^, 

jr(w) • n {a^p'^u-fi^ a~p~u-n^ pu-nui +pni, pu-nu2+pn2^ pHu-n) , (25) 

5(w) = {0,0, pgi,pg2,pg'u) . (26) 
The Jacobian matrix A(w) • n is deflned by 



In order to compute A(w) • n, one writes Eq. ( l25l ) for T^w) • n in terms of w 
and p: 

>(w) 'n = [wi ,W2 ■ ,ws hpni, 

V W1+W2 Wi + W2 Wi + W2 

W4 ■ +pn2, [w^+P) ■ . (28) 

Wi + W2 W1+W2 J 



The Jacobian matrix ( 1271 ) then has the following expression: 
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A(w) ■ n = 

-Un^ 

-UlUn + £-ni -UlUn + + umi + tima + l^ni 

-«2«n + -«2Wn + ^2^1 + ^n2 + ^2^2 + ^n2 ^n2 

.«n(fe-i^) u^{-t.-H) u^^^+Hm u^^^+Hn, «n(l + fe)/ 
where Un = u • n. 

Let us now compute the five derivatives dp/dwi. A systematic way of doing it 
is to introduce a set of five independent physical variables and here we shall 
take: 

(fl=a, (f2=P, ^3=T, (f4 = Ui, (f5=U2. (29) 

The expressions of the w[s in terms of the (fjS are algebraic and explicit. 
Hence the Jacobian matrix dwi/dipj can be easily computed. Since dipj/dwi 
is its inverse matrix, one finds easily with the help of a computer algebra 
program that 



{u^ + U2) + a p X , (30) 



dwi 2 

\m^ + «^) + «Vx+, (31) 



dw2 2 



where 



,.,,-.0, (33) 

7^ — 1 7^ — 1 
{cfr^ChHl^-l)T=^-^^, (34) 

r 

r-l^(7(a)-l) (35) 
In Eq. (l35l). we have introduced the speed of sound of the mixture c^, defined 

by 

J_ ^ (l + a)7+ (1 - «)7- _ J_ 

pc2 2p+(c+)2 ^ 2p-(c-)2 pa2 ' ^ ' 

with 

2 _ (l + a)p+(c+)^ (l-«)p-(c;)^ 
^ - 2(7+ - 1) + 2(7- - 1) J ^''^ 
Then one can show that the Jacobian matrix A(w) • n has three distinct 
eigenvalues: 

Xl=Un-Cs, A2,3,4 = ^n, A5 = li^ + C^, (38) 
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These three eigenvalues are real and there is a complete set of real valued 
eigenvectors. The expressions of these eigenvectors can be obtained by using 
a computer algebra program. 

Remark 2 //^+ = and tt' = 0, then d = ^ and = ^4t- 

Remark 3 The left hand side of l[SW is positive since pa^ is hounded from 
below by ) ^ (i-^^p Jc^ ) ^ j^/^^^ ^2 ^^^^ p^^y ^^^^ of the square 

of the enthalpy by analogy with the monofluid case. 

A plot of l/cs{a) is given in Fig. [1] (see the solid line). A remarkable property 
is that the speed of sound exhibits a minimum. If the energy equation was not 
taken into consideration, this minimum would not be present (see the dashed 
line in Fig. [1]). Indeed, the expression for the speed of sound for homogeneous 
two-fluid models with the additional assumption that the flow is isentropic is 



This expression can be found in Appendix C of [7j (see equations (C.20) and 



2.3 Evolution equations for the physical variables 

The system of conservation laws can be transformed into a set of 

evolution equations for the physical variables. Let us introduce the entropy 
function s{x^t) deflned by (compare with Eq. (fTOl) ) 



Cis 



{a-p^ + a^p-){ctY{c-Y 

\ a^p-{c-y + a-p^{cty ' 



(39) 



(C.21)). 



2p5 = (1 + a)p+s+ + (1 - a)p~s~. 



Proposition 1 Continuous solutions to (Cy ([^ satisfy 



Ut + u ' Vu H — = g 
P 



(40) 



Pt + u ' Vp + pclV ' u = 
at + u-Va + {l-a'^)SV -11=0 
St + u • Vs 



(41) 
(42) 
(43) 



where 



is given by ^3E)-li3^ and S is given by 



-2 p+p-(c+)2(c-)2 



1 pc^(7 7r+ — 7+7r ) 



(44) 
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qCI ^ ^ ^ ^ ^ ^ ^ ^ ^ 1 

-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 0.6 0.8 1 

OL 



(a) 




-1 -0.995 -0.99 



(b) 

Fig. 1. Speed of sound as a function of a. In order to make the curves more visible, 
the inverse of the speed of sound is plotted. The end points, represented by circles, 
correspond to pure liquid (left) and pure gas (right). The solid line represents the 
inverse of the speed of sound l/c^ given by equations (f36l) and (f37l) . The dashed line 
represents the inverse of the speed of sound l/c/^ (f39l) obtained without using the 
energy equation, (a) Full range; (b) Zoom near the pure liquid limit. The various 
parameters are: p+ = 1000 kg/m^, 7+ = 7, 7r+ = 2.1 x 10^ Pa, c+ = 1500 m/s and 
p- = 1.29 kg/m^ 7- = 1.4, TT- = Pa, c' = 300 m/s. 
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Remark 4 For pure fluids (a = ±1)^ Eq. is no longer relevant and 5 
is not needed. One can check that the speed of sound Cg is then equal to the 
expected speed of sound (of or c~ ) for pure fluids. 

The balance of entropy (H3l) comes from the balance 

{ps)t + V ■ {psu) = 0. (45) 

Adding together Eqs ^ and ([2]) leads to 

A + V • {pu) = 0. (46) 

Combining Eqs (USD and (SSI) leads to Eq. 

Remark 5 Subtracting Eq. ^ from Eq. ^ leads to 

{px)t + V . {pxu) = , with X = — . (47) 

P 

In the case of smooth solutions^ we obtain that 

+ ^ • Vx = , 
which is an alternative to Eq. [4^ - 

2.4 Pure fluid limit 

The two-fluid model described in the present paper is based on the volume 
fraction of liquid and gas. In some situations, this volume fraction can have 
sharp gradients. Consider for example a tanh-type distribution of a along the 
vertical axis with essentially pure gas at the top, pure liquid at the bottom 
and a middle layer where a goes rapidly from —1 to 1. One can even consider 
the limiting case where the transition is discontinuous. In this section we 
study this limit and we show that the two-fluid model degenerates into the 
classical water-wave equations. In other words one has an interface separating 
two pure fluids. So the well-known water-wave equations are a by-product of 
the two-fluid system under investigation. A similar type of limit in the case 
of a continuously stratifled incompressible fluid degenerating into a two-layer 
incompressible fluid was considered by James 

In the rest of this section, it is assumed that there are no shocks. Consider the 
3D case where a is either 1 or —1. More precisely let 

a := 1 — 27i{z — 7y(x, t)) , x = (xi, X2) , (48) 
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where Ti is the Heaviside step function, z the vertical coordinate and Xi,X2 
the horizontal coordinates. Physically this substitution means that we consider 
two pure fluids separated by an interface. It follows that 

a^a~ = , 1 — = . 



Substituting the expression (l48l) into the equation (l42l) gives 

r]t + Uh' ^hV = ^ , 

where Uh = (t^i, 1^2), = (^xi, ^0:2) w is the vertical velocity. 

This equation simply states that there is no mass flux across the interface. 
Incidentally this is no longer true in the case of shock waves. Integrating the 
conservation of momentum equation ([3j) inside a volume moving with the flow 
and enclosing the interface, and using the fact that there is no mass flux 
across the interface simply leads to the fact that there is no pressure jump 
across the interface. In other words, the pressure is continuous across the 
interface. Integrating the entropy equation inside the same volume enclosing 
the interface and using the fact there is no mass flux across the interface does 
not lead to any new information. 

One can now write Eqs ([2j)"(|ll) in each fluid by taking = 1, either in the 
conservative form 



{p% + V-ip^u^) = 0, (49) 
{p^u^)t + V ■ ip^u^®u^) + Vp^=p^g, (50) 
ip^s% + V-ip^s^u^) = 0, (51) 

(see Whitham [20] for example for the last equation) or in the more classical 
form 



+ {u^ . V)p^ + p^V'U^ = 0, (52) 
u^ + {u^.\/)u^ + ^ = g^ (53) 
+ • = . (54) 

In these two systems, the superscripts + and — are used for the heavy fluid 
(below the interface) and the light fluid (above the interface) respectively. 

The system of equations we derived is nothing else than the system of a 
discontinuous two-fluid system with an interface located at z = r]{x^ t). Along 
the interface, one has the kinematic and dynamic boundary conditions 



11 



r]t + u^ 'Vhr] = w^ , (55) 
p~=p'^. (56) 



This simple computation shows an important property of our modeh it auto- 
maticaUy degenerates into a discontinuous two-fluid system where two pure 
compressible phases are separated by an interface. This limit has interesting 
consequences. In particular, interfacial flows develop waves along the interface 
and these waves are usually dispersive. Therefore one can also expect disper- 
sive waves to exist in the two-fluid model. Since the emphasis of the present 
paper is the study of large-amplitude disturbances, the derivation of the dis- 
persion relation for the two-fluid model is left for future work. Note however 
that preliminary results can be found in [7j. Even the question of which rest 
state one must consider is not trivial. 



3 Simulations of aerated violent flows 



3.1 A finite-volume discretization of the model 



Here we describe the discretization of the model by a standard cell- 

centered flnite volume method. The computational domain C is trian- 
gulated into a set of control volumes: Q = UxerK. We start by integrating 
equation (l23l) on K: 

^ [ wdn+ y [ j^(w) . Hkl da= [ 5(w) dQ , (57) 

dt JK L^{K) -^^^^ 

where Hkl denotes the unit normal vector on K H L pointing into L and 
X{K) = {LeT : area(X n L) 7^ 0} . Then, setting 



WK{t) := ] [ w(x,t) dn , 

vol A JK 



we approximate ( l57l ) by 

dwK area(Lni^)^. ^ . . 

where the numerical flux 

$(wk, w^; ukl) ~ ^ J.. I ^(w) ■ ukl da 

area(L fl K) Jkhl 
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is explicitly computed by the FVCF formula of Ghidaglia et al. [H] : 



$(v, w; n) = ^ sgn(A„(/i(v, w)))^^ . 

(59) 

Here v and w are dummy variables. The Jacobian matrix Anda) is defined in 
( l27l ). /x(v,w) is an arbitrary mean between v and w (for example /x(v,w) = 
(l/2)(v+w)) and sgn(M) is the matrix whose eigenvectors are those of M but 
whose eigenvalues are the signs of that of M. As explained in [9j, this method 
is able to model discontinuities such as shock waves and sharp interfaces. 

So far we have not discussed the case where a control volume K meets the 
boundary of Q. Here we shall only consider the case where this boundary is a 
wall and from the numerical point of view, we only need to find the normal 
flux T • n. Since u{x^ t) • n = for x G , we have 

i^'^)\xedn = (0,0,^5^,0), := , 

and following Ghidaglia and Pascal [10], we can take Pb = P + pUnCs^ where 
the right-hand side is evaluated in the control volume K. 

Remark 1 In order to turn l[5W into a numerical algorithm^ we must at 
least perform time discretization and give an expression for /x(v,w). Since 
this matter is standard^ we do not give the details here but instead refer to 
Dutykh f^J. Let us also notice that formula (E^j leads to a first- order scheme 
but in fact we use a MUSCL technique to achieve higher accuracy in space 

m- 

3.2 Numerical results 

In order to check the accuracy of our second-order scheme on smooth solutions 
and its robustness against discontinuous solutions, we have performed the 
classical test cases for which we refer to [8] . The most famous test case is that 
of Sod's shock tube. We report here on some of the situations which have 
motivated this study. 

3.2.1 Thermodynamics constants 

The constants Cy can be calculated after simple algebraic manipulations of 
equations ([71), ([8]) and matching with experimental values at normal condi- 
tions: 

r- = 

^" (7--iWro' 
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parameter value 



Po 



10^ Pa 



Pq 10^ kg/m^ 

Po 1.29 kg/rn^ 

To 300 K 

7" 1.4 



7r+ 2.1 X 10^ Pa 
C,7 646.0 



Table 1 

Values of the parameters for an air/water mixture under normal conditions. 

^+ ^ 7+Po + 7r+ 
^ - (7+ - l)7+P^ro ■ 
For example, for an air/water mixture under normal conditions we have the 
values given in Table [H 

The sound velocities in each phase are given by the following formulas: 

ic:r = (cif = '-^S^. (60) 



In the two test cases described below, we use a very high value for the accel- 
eration due to gravity: g = 100 ms~^. The only motivation is to accelerate the 
dynamics. All results are presented with physical dimensions. For example, 
the 1x1 box used for the computations corresponds to a 1 m by 1 m box. 



3.2.2 Falling water column 

The geometry and initial condition for this test case are shown on Fig. [21 
Initially the velocity field is taken to be zero. At time t = 0, the volume 
fraction of gas is 0.9 (white area) while the volume fraction of water is 0.9 
(dark area) . The values of the other parameters are given in Table [H The 
mesh used in this computation contained about 108000 control volumes (in 
this case they were triangles). The results of this simulation are presented on 
Figures [HHSl Fig. [9] shows the maximal pressure on the right wall as a function 
of time: 

t I — > max p(x, y, t). 
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We performed another computation for a mixture with a"*" = 0.05, a~ = 0.95. 
The pressure is recorded as weU and plotted in Fig. [TUl One can see that the 
peak value is higher and the impact is more localized in time. 



1 



0.9 



«+ = 0.9 
a- = 0.1 







0.3 



tt+ = 0.1 
a- = 0.9 



Hi 



0.05 



0.65 0.7 



Fig. 2. Falling water column test case. Geometry and initial condition. All the values 
for are at time t = 0. 



3.2.3 Water drop test case 



The geometry and initial condition for this test case are shown on Fig. [TTl 
Initially the velocity field is taken to be zero. The values of the other parame- 




Fully compr3S5itale homogeneous t'vo phase soh/er 
Mixture density at t = 0.005 



Fully compr3S5itale homogeneous tw pf-ase soh/er 
Mixture density at t = 0.060 




(a) t = 0.005 s 



(b) t = 0.06 s 



Fig. 3. Falling water column test case. Initial condition and the beginning of the 
column collapse. 
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Fully c3mpr3S5ibte homogeneous Ivo phase soh/er 
Mixture density at t = 0.100 
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Fully c3mpr3S5ible homogeneous t'vo phase soh/er 
Mixture density at t = 0.125 




(a) t = 0.1 s 



(b) t = 0.125 s 



Fig. 4. Falling water column test case. Splash formation due to the interaction with 
the step. 



Fully c5mpr3S5itale homogeneous [vo phase soh/er 
Mixture density at t = 0.150 




Fully compr3S5itale homogeneous tw phase soh/er 
Mixture density at t = 0.175 




(a) t = 0.15 s 



(b) t = 0.175 s 



Fig. 5. Falling water column test case. Water hits the wall. 



Fully compr3S5itale homogeneous [vo phase soh/er 
Mixture density at t = 0.200 




Fully compr3S5itale homogeneous tw phase soh/er 
Mixture density at t = 0.225 




(a) t = 0.2 s (b) t = 0.225 s 

Fig. 6. Same as Fig. Oat later times. 



ters are given in Table [H The mesh used in this computation contained about 
92000 control volumes (again they were triangles). The results of this simula- 
tion are presented in Figures [T214T81 In Fig. [19] we plot the maximal pressure 
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Fully csmprsssibte homogeneous t'w p|-ase soh/er 
Mixture density at t = 0.300 




Fully comprsssible homogeneous t'w phase soh/er 
Mixture density at t = 0.400 



C.2 0.4 o.s 0.8 




(a) t = 0.3 s (b) t = 0.4 s 

Fig. 7. Falling water column test case. The splash is climbing the wall. 



Fully comprsssitale homogeneous t'vo phase soh/er 
Mixture density at t = 0.500 




Fully comprsssitale homogeneous Ivo phase soh/er 
Mixture density at t = 0.675 




(a) t = 0.5 s (b) t = 0.675 s 

Fig. 8. Falling water column test case. Turbulent mixing process. 
Maximal pressure on the right wall 




0.3 0.4 
t, time 



Fig. 9. Maximal pressure on the right wall as a function of time. Case of a heavy 
gas. 



17 



Maximal pressure on the right wall 



2.2 



E 
Q. 



0.8 









































































J 

























0.1 0.2 0.3 0.4 0.5 0.6 0.7 

t, time 



Fig. 10. Maximal pressure on the right wall as a function of time. Case of a light 
gas. 

on the bottom as a function of time: 



t I — > max p(x, y, t). 

(x,y)G[0,l]xO 



The pressure exerted on the bottom reaches 2.5po due to the drop impact at 
t ^ 0.16 s. 



4 Conclusions 



In this article we have presented a simple mathematical model for simulating 
water wave impacts. Associated to this model, which avoids the costly capture 
of free surfaces, we have built a numerical solver which is: (i) second-order 
accurate on smooth solutions, (ii) stable even for solutions with very strong 
gradients (and solutions with shocks) and (iii) locally exactly conservative 
with respect to the mass of each fluid, momentum and total energy. This last 
property, (iii), which is certainly the most desirable from the physical point of 
view, is an immediate byproduct of our cell-centered flnite volume method. 

We have shown here the good behavior of this framework on simple test cases 
and we are presently working on quantitative comparisons in the context of 
real applications. 
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Fig. 11. Geometry and initial condition for water drop test case. All the values for 
azb are at time t = 0. 



Mixtue densty at t =0.005 



MiKttre densty at t =0.075 
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(a) t = 0.005 s (b) t = 0.075 s 

Fig. 12. Water drop test case. Initial configuration and the beginning of the fall. 
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Mixture d«ns ty at t = .1 03 



MiKtLre dens ty at t = .1 25 
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(a) t = 0.1 s (b) t = 0.125 s 

Fig. 13. Water drop test case. Drop approaching the bottom of the container. 



MiKtue densty at t = 0.1 35 



MiKttre densty at t = 0.1 53 



0.2 0.4 O.e 



0.2 0.4 O.e o.a 



(a) t = 0.135 s (b) t = 0.15 s 

Fig. 14. Water drop test case. Drop/bottom compressible interaction. 
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MiKtue densty at t =0.203 
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(b) t = 0.2 s 



Fig. 15. Water drop test case. Formation of vertical jets. 
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Fig. 16. Water drop test case. Crossing of side jets. 
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